home *** CD-ROM | disk | FTP | other *** search
/ Aminet 35 / Aminet 35 (2000)(Schatztruhe)[!][Feb 2000].iso / Aminet / gfx / misc / gnuplot-src.lha / gnuplot-3.7.1src / gnuplot-3.7.1.lha / gnuplot-3.7.1 / plot2d.c < prev    next >
Encoding:
C/C++ Source or Header  |  1999-10-11  |  43.7 KB  |  1,421 lines
  1. #ifndef lint_AXIS])  = yp;
  2.     Sid = "$Id: plot2d.c,v 1.16.2.4 1999/10/11 12:18:56 lhecking Exp $";
  3. #endif
  4.  
  5. /* GNUPLOT - plot2d.c */
  6.  
  7. /*[
  8.  * Copyright 1986 - 1993, 1998   Thomas Williams, Colin Kelley
  9.  *
  10.  * Permission to use, copy, and distribute this software and its
  11.  * documentation for any purpose with or without fee is hereby granted,
  12.  * provided that the above copyright notice appear in all copies and
  13.  * that both that copyright notice and this permission notice appear
  14.  * in supporting documentation.
  15.  *
  16.  * Permission to modify the software is granted, but not the right to
  17.  * distribute the complete modified source code.  Modifications are to
  18.  * be distributed as patches to the released version.  Permission to
  19.  * distribute binaries produced by compiling modified sources is granted,
  20.  * provided you
  21.  *   1. distribute the corresponding source modifications from the
  22.  *    released version in the form of a patch file along with the binaries,
  23.  *   2. add special version identification to distinguish your version
  24.  *    in addition to the base release version number,
  25.  *   3. provide your name and address as the primary contact for the
  26.  *    support of your modified version, and
  27.  *   4. retain our contact information in regard to use of the base
  28.  *    software.
  29.  * Permission to distribute the released version of the source code along
  30.  * with corresponding source modifications in the form of a patch file is
  31.  * granted with same provisions 2 through 4 for binary distributions.
  32.  *
  33.  * This software is provided "as is" without express or implied warranty
  34.  * to the extent permitted by applicable law.
  35. ]*/
  36.  
  37. #include "plot.h"
  38. #include "setshow.h"
  39. #include "fit.h"
  40. #include "binary.h"
  41.  
  42. #ifndef _Windows
  43. # include "help.h"
  44. #endif
  45.  
  46. #ifndef STDOUT
  47. #define STDOUT 1
  48. #endif
  49.  
  50. /* static prototypes */
  51.  
  52. void plotrequest __PROTO((void));
  53. void plot3drequest __PROTO((void));
  54. void define __PROTO((void));
  55. static void get_data __PROTO((struct curve_points * this_plot));
  56. static void store2d_point __PROTO((struct curve_points * this_plot, int i, double x, double y, double xlow, double xhigh, double ylow, double yhigh, double width));
  57. static void print_table __PROTO((struct curve_points * first_plot, int plot_num));
  58. static void eval_plots __PROTO((void));
  59. static void parametric_fixup __PROTO((struct curve_points * start_plot, int *plot_num));
  60.  
  61.  
  62. /* the curves/surfaces of the plot */
  63. struct curve_points *first_plot = NULL;
  64. static struct udft_entry plot_func;
  65. extern struct udft_entry *dummy_func;
  66.  
  67. /* jev -- for passing data thru user-defined function */
  68. /* Needed by datafile.c */
  69. struct udft_entry ydata_func;
  70.  
  71. extern int datatype[];
  72. extern char timefmt[];
  73.  
  74. extern TBOOLEAN is_3d_plot;
  75. extern int plot_token;
  76.  
  77. /* in order to support multiple axes, and to
  78.  * simplify ranging in parametric plots, we use
  79.  * arrays to store some things.
  80.  * Elements are z = 0, y1 = 1, x1 = 2, [z2 =4 ], y2 = 5, x2 = 6
  81.  * these are given symbolic names in plot.h
  82.  */
  83.  
  84. /* Were declared in command.c */
  85. double min_array[AXIS_ARRAY_SIZE], max_array[AXIS_ARRAY_SIZE];
  86. int auto_array[AXIS_ARRAY_SIZE];
  87. TBOOLEAN log_array[AXIS_ARRAY_SIZE];
  88. double base_array[AXIS_ARRAY_SIZE];
  89. double log_base_array[AXIS_ARRAY_SIZE];
  90.  
  91. /* Deleted from setshow.h and renamed */
  92. extern FILE *gpoutfile;
  93.  
  94. /* if user specifies [10:-10] we use [-10:10] internally, and swap at end */
  95. int reverse_range[AXIS_ARRAY_SIZE];
  96.  
  97. /* info from datafile module */
  98. extern int df_datum;
  99. extern int df_line_number;
  100. extern int df_no_use_specs;
  101. extern int df_eof;
  102. extern int df_timecol[];
  103. extern TBOOLEAN df_binary;
  104.  
  105. #define Inc_c_token if (++c_token >= num_tokens)    \
  106.                         int_error ("Syntax error", c_token);
  107.  
  108.  
  109. /*
  110.  * IMHO, code is getting too cluttered with repeated chunks of
  111.  * code. Some macros to simplify, I hope.
  112.  *
  113.  * do { } while(0) is comp.lang.c recommendation for complex macros
  114.  * also means that break can be specified as an action, and it will
  115.  * 
  116.  */
  117.  
  118. /*  copy scalar data to arrays
  119.  * optimiser should optimise infinite away
  120.  * dont know we have to support ranges [10:-10] - lets reverse
  121.  * it for now, then fix it at the end.
  122.  */
  123. #define INIT_ARRAYS(axis, min, max, auto, is_log, base, log_base, infinite) \
  124. do{auto_array[axis] = auto; \
  125.    min_array[axis] = (infinite && (auto&1)) ? VERYLARGE : min; \
  126.    max_array[axis] = (infinite && (auto&2)) ? -VERYLARGE : max; \
  127.    log_array[axis] = is_log; base_array[axis] = base; log_base_array[axis] = log_base;\
  128. }while(0)
  129. /* handle reversed ranges */
  130. #define CHECK_REVERSE(axis) \
  131. do{\
  132.  if (auto_array[axis] == 0 && max_array[axis] < min_array[axis]) {\
  133.   double temp = min_array[axis]; min_array[axis] = max_array[axis]; max_array[axis] = temp;\
  134.   reverse_range[axis] = 1; \
  135.  } else reverse_range[axis] = (range_flags[axis]&RANGE_REVERSE); \
  136. }while(0)
  137.  
  138.  
  139. /* get optional [min:max] */
  140. #define LOAD_RANGE(axis) \
  141. do {\
  142.  if (equals(c_token, "[")) { \
  143.   c_token++; \
  144.   auto_array[axis] = load_range(axis,&min_array[axis], &max_array[axis], auto_array[axis]);\
  145.   if (!equals(c_token, "]"))\
  146.    int_error("']' expected", c_token);\
  147.   c_token++;\
  148.  }\
  149. } while (0)
  150.  
  151.  
  152. /* store VALUE or log(VALUE) in STORE, set TYPE as appropriate
  153.  * Do OUT_ACTION or UNDEF_ACTION as appropriate
  154.  * adjust range provided type is INRANGE (ie dont adjust y if x is outrange
  155.  * VALUE must not be same as STORE
  156.  */
  157.  
  158. #define STORE_WITH_LOG_AND_FIXUP_RANGE(STORE, VALUE, TYPE, AXIS, OUT_ACTION, UNDEF_ACTION)\
  159. do { if (log_array[AXIS]) { if (VALUE<0.0) {TYPE = UNDEFINED; UNDEF_ACTION; break;} \
  160.               else if (VALUE == 0.0){STORE = -VERYLARGE; TYPE = OUTRANGE; OUT_ACTION; break;} \
  161.               else { STORE = log(VALUE)/log_base_array[AXIS]; } \
  162.      } else STORE = VALUE; \
  163.      if (TYPE != INRANGE) break;  /* dont set y range if x is outrange, for example */ \
  164.      if ( VALUE<min_array[AXIS] ) { \
  165.       if (auto_array[AXIS] & 1) min_array[AXIS] = VALUE; else { TYPE = OUTRANGE; OUT_ACTION; break; }  \
  166.      } \
  167.      if ( VALUE>max_array[AXIS] ) { \
  168.       if (auto_array[AXIS] & 2) max_array[AXIS] = VALUE; else { TYPE = OUTRANGE; OUT_ACTION; }   \
  169.      } \
  170. } while(0)
  171.  
  172. /* use this instead empty macro arguments to work around NeXT cpp bug */
  173. /* if this fails on any system, we might use ((void)0) */
  174. #define NOOP            /* */
  175.  
  176. /* check range and take logs of min and max if logscale
  177.  * this also restores min and max for ranges like [10:-10]
  178.  */
  179. #ifdef HAVE_STRINGIZE
  180. # define LOG_MSG(x) #x " range must be greater than 0 for log scale!"
  181. #else
  182. # define LOG_MSG(x) "x range must be greater than 0 for log scale!"
  183. #endif
  184.  
  185. #define FIXUP_RANGE_FOR_LOG(AXIS, WHICH) \
  186. do { if (reverse_range[AXIS]) { \
  187.       double temp = min_array[AXIS]; \
  188.       min_array[AXIS] = max_array[AXIS]; \
  189.       max_array[AXIS] = temp; \
  190.      }\
  191.      if (log_array[AXIS]) { \
  192.       if (min_array[AXIS] <= 0.0 || max_array[AXIS] <= 0.0) \
  193.        int_error(LOG_MSG(WHICH), NO_CARET); \
  194.       min_array[AXIS] = log(min_array[AXIS])/log_base_array[AXIS]; \
  195.       max_array[AXIS] = log(max_array[AXIS])/log_base_array[AXIS];  \
  196.    } \
  197. } while(0)
  198.  
  199.  
  200.  
  201. void plotrequest()
  202. /*
  203.  * In the parametric case we can say plot [a= -4:4] [-2:2] [-1:1] sin(a),a**2
  204.  * while in the non-parametric case we would say only plot [b= -2:2] [-1:1]
  205.  * sin(b)
  206.  */
  207. {
  208.     int dummy_token = -1;
  209.  
  210.     if (!term)            /* unknown */
  211.     int_error("use 'set term' to set terminal type first", c_token);
  212.  
  213.     is_3d_plot = FALSE;
  214.  
  215.     if (parametric && strcmp(dummy_var[0], "u") == 0)
  216.     strcpy(dummy_var[0], "t");
  217.  
  218.     /* initialise the arrays from the 'set' scalars */
  219.  
  220.     INIT_ARRAYS(FIRST_X_AXIS, xmin, xmax, autoscale_x, is_log_x, base_log_x, log_base_log_x, 0);
  221.     INIT_ARRAYS(FIRST_Y_AXIS, ymin, ymax, autoscale_y, is_log_y, base_log_y, log_base_log_y, 1);
  222.     INIT_ARRAYS(SECOND_X_AXIS, x2min, x2max, autoscale_x2, is_log_x2, base_log_x2, log_base_log_x2, 0);
  223.     INIT_ARRAYS(SECOND_Y_AXIS, y2min, y2max, autoscale_y2, is_log_y2, base_log_y2, log_base_log_y2, 1);
  224.  
  225.     min_array[T_AXIS] = tmin;
  226.     max_array[T_AXIS] = tmax;
  227.  
  228.     if (equals(c_token, "[")) {
  229.     c_token++;
  230.     if (isletter(c_token)) {
  231.         if (equals(c_token + 1, "=")) {
  232.         dummy_token = c_token;
  233.         c_token += 2;
  234.         } else {
  235.         /* oops; probably an expression with a variable. */
  236.         /* Parse it as an xmin expression. */
  237.         /* used to be: int_error("'=' expected",c_token); */
  238.         }
  239.     } {
  240.         int axis = (parametric || polar) ? T_AXIS : FIRST_X_AXIS;
  241.  
  242.  
  243.         auto_array[axis] = load_range(axis, &min_array[axis], &max_array[axis], auto_array[axis]);
  244.         if (!equals(c_token, "]"))
  245.         int_error("']' expected", c_token);
  246.         c_token++;
  247.     }            /* end of scope of 'axis' */
  248.     }                /* first '[' */
  249.     if (parametric || polar)    /* set optional x ranges */
  250.     LOAD_RANGE(FIRST_X_AXIS);
  251.     else {
  252.     /* order of t doesn't matter, but x does */
  253.     CHECK_REVERSE(FIRST_X_AXIS);
  254.     }
  255.  
  256.     LOAD_RANGE(FIRST_Y_AXIS);
  257.     CHECK_REVERSE(FIRST_Y_AXIS);
  258.     LOAD_RANGE(SECOND_X_AXIS);
  259.     CHECK_REVERSE(SECOND_X_AXIS);
  260.     LOAD_RANGE(SECOND_Y_AXIS);
  261.     CHECK_REVERSE(SECOND_Y_AXIS);
  262.  
  263.     /* use the default dummy variable unless changed */
  264.     if (dummy_token >= 0)
  265.     copy_str(c_dummy_var[0], dummy_token, MAX_ID_LEN);
  266.     else
  267.     (void) strcpy(c_dummy_var[0], dummy_var[0]);
  268.  
  269.     eval_plots();
  270. }
  271.  
  272. /* Use up to 7 columns in data file at once -- originally it was 5 */
  273. #define NCOL 7
  274.  
  275.  
  276. /* A quick note about boxes style. For boxwidth auto, we cannot
  277.  * calculate widths yet, since it may be sorted, etc. But if
  278.  * width is set, we must do it now, before logs of xmin/xmax
  279.  * are taken.
  280.  * We store -1 in point->z as a marker to mean width needs to be
  281.  * calculated, or 0 to mean that xmin/xmax are set correctly
  282.  */
  283.  
  284.  
  285. static void get_data(this_plot)
  286. struct curve_points *this_plot;
  287. /* this_plot->token is after datafile spec, for error reporting
  288.  * it will later be moved passed title/with/linetype/pointtype
  289.  */
  290. {
  291.     int i /* num. points ! */ , j;
  292.     int max_cols, min_cols;    /* allowed range of column numbers */
  293.     double v[NCOL];
  294.     int storetoken = this_plot->token;
  295.  
  296.     /* eval_plots has already opened file */
  297.  
  298.     switch (this_plot->plot_style) {    /* set maximum columns to scan */
  299.     case XYERRORBARS:
  300.     case BOXXYERROR:
  301.     max_cols = 7;
  302.     min_cols = 4;
  303.     break;
  304.  
  305.     case FINANCEBARS:
  306.     case CANDLESTICKS:
  307.     min_cols = max_cols = 5;
  308.     break;
  309.  
  310.     case BOXERROR:
  311.     max_cols = 5;
  312.     min_cols = 4;
  313.     break;
  314.  
  315.     case VECTOR:
  316.     min_cols = max_cols = 4;
  317.     break;
  318.  
  319.     case XERRORBARS:
  320.     case YERRORBARS:
  321.     max_cols = 4;
  322.     min_cols = 3;
  323.     break;
  324.  
  325.     case BOXES:
  326.     min_cols = 2;
  327.     max_cols = 4;
  328.  
  329.     default:
  330.     min_cols = 1;
  331.     max_cols = 2;
  332.     }
  333.  
  334.     if (this_plot->plot_smooth == ACSPLINES)
  335.     max_cols = 3;
  336.  
  337.     if (df_no_use_specs > max_cols)
  338.     int_error("Too many using specs for this style", NO_CARET);
  339.  
  340.     if (df_no_use_specs >0 && df_no_use_specs < min_cols)
  341.     int_error("Not enough columns for this style", NO_CARET);
  342.  
  343.     i = 0;
  344.     while ((j = df_readline(v, max_cols)) != DF_EOF) {
  345.     /* j <= max_cols */
  346.  
  347.     if (i >= this_plot->p_max) {
  348.         /*
  349.          * overflow about to occur. Extend size of points[] array. We
  350.          * either double the size, or add 1000 points, whichever is a
  351.          * smaller increment. Note i = p_max.
  352.          */
  353.         cp_extend(this_plot, i + (i < 1000 ? i : 1000));
  354.     }
  355.     /* Limitation: No xerrorbars with boxes */
  356.     switch (j) {
  357.     default:
  358.         {
  359.         char message[80];
  360.         sprintf(message, "internal error : df_readline returned %d : datafile line %d", j, df_line_number);
  361.         df_close();
  362.         int_error(message, c_token);
  363.         }
  364.     case DF_UNDEFINED:
  365.         /* bad result from extended using expression */
  366.         this_plot->points[i].type = UNDEFINED;
  367.         i++;
  368.         continue;
  369.  
  370.     case DF_FIRST_BLANK:
  371.         /* break in data, make next point undefined */
  372.         this_plot->points[i].type = UNDEFINED;
  373.         i++;
  374.         continue;
  375.  
  376.     case DF_SECOND_BLANK:
  377.         /* second blank line. We dont do anything
  378.          * (we did everything when we got FIRST one)
  379.          */
  380.         continue;
  381.  
  382.     case 0:        /* not blank line, but df_readline couldn't parse it */
  383.         {
  384.         char message[80];
  385.         sprintf(message, "Bad data on line %d", df_line_number);
  386.         df_close();
  387.         int_error(message, this_plot->token);
  388.         }
  389.  
  390.     case 1:
  391.         {            /* only one number */
  392.         /* x is index, assign number to y */
  393.         v[1] = v[0];
  394.         v[0] = df_datum;
  395.         /* nobreak */
  396.         }
  397.  
  398.     case 2:
  399.         /* x, y */
  400.         /* ylow and yhigh are same as y */
  401.  
  402.         if (this_plot->plot_style == BOXES && boxwidth > 0) {
  403.         /* calc width now */
  404.         store2d_point(this_plot, i++, v[0], v[1], v[0] - boxwidth / 2, v[0] + boxwidth / 2, v[1], v[1], 0.0);
  405.         } else {
  406.         /* xlow and xhigh are same as x */
  407.         store2d_point(this_plot, i++, v[0], v[1], v[0], v[0], v[1], v[1],
  408.                   -1.0);    /* auto width if boxes, else ignored */
  409.         }
  410.         break;
  411.  
  412.  
  413.     case 3:
  414.         /* x, y, ydelta OR x, y, xdelta OR x, y, width */
  415.         if (this_plot->plot_smooth == ACSPLINES)
  416.         store2d_point(this_plot, i++, v[0], v[1], v[0], v[0], v[1], v[1], v[2]);
  417.         else
  418.         switch (this_plot->plot_style) {
  419.         default:
  420.             int_warn("This plot style not work with 3 cols. Setting to yerrorbars", storetoken);
  421.             this_plot->plot_style = YERRORBARS;
  422.             /* fall through */
  423.  
  424.         case YERRORBARS:
  425.         case BOXERROR:    /* x, y, dy */
  426.             store2d_point(this_plot, i++, v[0], v[1], v[0], v[0], v[1] - v[2], v[1] + v[2],
  427.                   -1.0);    /* auto width if boxes, else ignored */
  428.             break;
  429.  
  430.         case XERRORBARS:
  431.             store2d_point(this_plot, i++, v[0], v[1], v[0] - v[2], v[0] + v[2], v[1], v[1], 0.0);
  432.             break;
  433.  
  434.         case BOXES:
  435.             /* calculate xmin and xmax here, so that logs are taken if
  436.              * if necessary
  437.              */
  438.             store2d_point(this_plot, i++, v[0], v[1], v[0] - v[2] / 2, v[0] + v[2] / 2, v[1], v[1], 0.0);
  439.             break;
  440.  
  441.         }        /*inner switch */
  442.  
  443.         break;
  444.  
  445.  
  446.  
  447.     case 4:
  448.         /* x, y, ylow, yhigh OR
  449.          * x, y, xlow, xhigh OR
  450.          * x, y, xdelta, ydelta OR
  451.          * x, y, ydelta, width
  452.          */
  453.  
  454.         switch (this_plot->plot_style) {
  455.         default:
  456.         int_warn("This plot style does not work with 4 cols. Setting to yerrorbars",
  457.              storetoken);
  458.         this_plot->plot_style = YERRORBARS;
  459.         /* fall through */
  460.  
  461.         case YERRORBARS:
  462.         store2d_point(this_plot, i++, v[0], v[1], v[0], v[0], v[2], v[3], -1.0);
  463.         break;
  464.  
  465.         case BOXXYERROR:    /* x, y, dx, dy */
  466.         case XYERRORBARS:
  467.         store2d_point(this_plot, i++, v[0], v[1], v[0] - v[2], v[0] + v[2], v[1] - v[3], v[1] + v[3], 0.0);
  468.         break;
  469.  
  470.  
  471.         case BOXES:    /* x, y, xmin, xmax */
  472.         store2d_point(this_plot, i++, v[0], v[1], v[2], v[3], v[1], v[1], 0.0);
  473.         break;
  474.  
  475.         case XERRORBARS:
  476.         store2d_point(this_plot, i++, v[0], v[1], v[2], v[3], v[1], v[1], 0.0);
  477.         break;
  478.  
  479.         case BOXERROR:
  480.         /* x,y, xleft, xright */
  481.         store2d_point(this_plot, i++, v[0], v[1], v[0], v[0], v[1] - v[2], v[1] + v[2], 0.0);
  482.         break;
  483.  
  484.         case VECTOR:
  485.         /* x,y,dx,dy */
  486.         store2d_point(this_plot, i++, v[0], v[1], v[0], v[0] + v[2], v[1], v[1] + v[3], -1.0);
  487.         break;
  488.         }            /*inner switch */
  489.  
  490.         break;
  491.  
  492.  
  493.     case 5:
  494.         {            /* x, y, ylow, yhigh, width  or  x open low high close */
  495.         switch (this_plot->plot_style) {
  496.         default:
  497.             int_warn("Five col. plot style must be boxerrorbars, financebars or candlesticks. Setting to boxerrorbars", storetoken);
  498.             this_plot->plot_style = BOXERROR;
  499.             /*fall through */
  500.  
  501.         case BOXERROR:    /* x, y, ylow, yhigh, width */
  502.             store2d_point(this_plot, i++, v[0], v[1], v[0] - v[4] / 2, v[0] + v[4] / 2, v[2], v[3], 0.0);
  503.             break;
  504.  
  505.         case FINANCEBARS:
  506.         case CANDLESTICKS:
  507.             store2d_point(this_plot, i++, v[0], v[1], v[0], v[0], v[2], v[3], v[4]);
  508.             break;
  509.         }
  510.         break;
  511.         }
  512.  
  513.     case 7:
  514.         /* same as six columns. Width ignored */
  515.         /* eh ? - fall through */
  516.     case 6:
  517.         /* x, y, xlow, xhigh, ylow, yhigh */
  518.         switch (this_plot->plot_style) {
  519.         default:
  520.         int_warn("This plot style not work with 6 cols. Setting to xyerrorbars", storetoken);
  521.         this_plot->plot_style = XYERRORBARS;
  522.         /*fall through */
  523.         case XYERRORBARS:
  524.         case BOXXYERROR:
  525.         store2d_point(this_plot, i++, v[0], v[1], v[2], v[3], v[4], v[5], 0.0);
  526.         break;
  527.         }
  528.  
  529.     }            /*switch */
  530.  
  531.     }                /*while */
  532.  
  533.     this_plot->p_count = i;
  534.     cp_extend(this_plot, i);    /* shrink to fit */
  535.  
  536.     df_close();
  537. }
  538.  
  539. /* called by get_data for each point */
  540. static void store2d_point(this_plot, i, x, y, xlow, xhigh, ylow, yhigh, width)
  541. struct curve_points *this_plot;
  542. int i;                /* point number */
  543. double x, y;
  544. double ylow, yhigh;
  545. double xlow, xhigh;
  546. double width;            /* -1 means autocalc, 0 means use xmin/xmax */
  547. {
  548.     struct coordinate GPHUGE *cp = &(this_plot->points[i]);
  549.     int dummy_type = INRANGE;    /* sometimes we dont care about outranging */
  550.  
  551.  
  552.     /* jev -- pass data values thru user-defined function */
  553.     /* div -- y is dummy variable 2 - copy value there */
  554.     if (ydata_func.at) {
  555.     struct value val;
  556.  
  557.     (void) Gcomplex(&ydata_func.dummy_values[0], y, 0.0);
  558.     ydata_func.dummy_values[2] = ydata_func.dummy_values[0];
  559.     evaluate_at(ydata_func.at, &val);
  560.     y = undefined ? 0.0 : real(&val);
  561.  
  562.     (void) Gcomplex(&ydata_func.dummy_values[0], ylow, 0.0);
  563.     ydata_func.dummy_values[2] = ydata_func.dummy_values[0];
  564.     evaluate_at(ydata_func.at, &val);
  565.     ylow = undefined ? 0 : real(&val);
  566.  
  567.     (void) Gcomplex(&ydata_func.dummy_values[0], yhigh, 0.0);
  568.     ydata_func.dummy_values[2] = ydata_func.dummy_values[0];
  569.     evaluate_at(ydata_func.at, &val);
  570.     yhigh = undefined ? 0 : real(&val);
  571.     }
  572.     dummy_type = cp->type = INRANGE;
  573.  
  574.     if (polar) {
  575.     double newx, newy;
  576.     if (!(autoscale_r & 2) && y > rmax) {
  577.         cp->type = OUTRANGE;
  578.     }
  579.     if (!(autoscale_r & 1)) {
  580.         /* we store internally as if plotting r(t)-rmin */
  581.         y -= rmin;
  582.     }
  583.     newx = y * cos(x * ang2rad);
  584.     newy = y * sin(x * ang2rad);
  585. #if 0 /* HBB 981118: added polar errorbars */
  586.     /* only lines and points supported with polar */
  587.     y = ylow = yhigh = newy;
  588.     x = xlow = xhigh = newx;
  589. #else 
  590.     y = newy;
  591.     x = newx;
  592.  
  593.     if (!(autoscale_r & 2) && yhigh > rmax) {
  594.         cp->type = OUTRANGE;
  595.     }
  596.     if (!(autoscale_r & 1)) {
  597.         /* we store internally as if plotting r(t)-rmin */
  598.         yhigh -= rmin;
  599.     }
  600.     newx = yhigh * cos(xhigh * ang2rad);
  601.     newy = yhigh * sin(xhigh * ang2rad);
  602.     yhigh = newy;
  603.     xhigh = newx;
  604.  
  605.     if (!(autoscale_r & 2) && ylow > rmax) {
  606.         cp->type = OUTRANGE;
  607.     }
  608.     if (!(autoscale_r & 1)) {
  609.         /* we store internally as if plotting r(t)-rmin */
  610.         ylow -= rmin;
  611.     }
  612.     newx = ylow * cos(xlow * ang2rad);
  613.     newy = ylow * sin(xlow * ang2rad);
  614.     ylow = newy;
  615.     xlow = newx;
  616. #endif
  617.     }
  618.     /* return immediately if x or y are undefined
  619.      * we dont care if outrange for high/low.
  620.      * BUT if high/low undefined (ie log( < 0 ), no number is stored,
  621.      * but graphics.c doesn't know.
  622.      * explicitly store -VERYLARGE;
  623.      */
  624.     STORE_WITH_LOG_AND_FIXUP_RANGE(cp->x, x, cp->type, this_plot->x_axis, NOOP, return);
  625.     STORE_WITH_LOG_AND_FIXUP_RANGE(cp->xlow, xlow, dummy_type, this_plot->x_axis, NOOP, cp->xlow = -VERYLARGE);
  626.     STORE_WITH_LOG_AND_FIXUP_RANGE(cp->xhigh, xhigh, dummy_type, this_plot->x_axis, NOOP, cp->xhigh = -VERYLARGE);
  627.     STORE_WITH_LOG_AND_FIXUP_RANGE(cp->y, y, cp->type, this_plot->y_axis, NOOP, return);
  628.     STORE_WITH_LOG_AND_FIXUP_RANGE(cp->ylow, ylow, dummy_type, this_plot->y_axis, NOOP, cp->ylow = -VERYLARGE);
  629.     STORE_WITH_LOG_AND_FIXUP_RANGE(cp->yhigh, yhigh, dummy_type, this_plot->y_axis, NOOP, cp->yhigh = -VERYLARGE);
  630.     cp->z = width;
  631. }                /* store2d_point */
  632.  
  633.  
  634.  
  635. /*
  636.  * print_points: a debugging routine to print out the points of a curve, and
  637.  * the curve structure. If curve<0, then we print the list of curves.
  638.  */
  639.  
  640. #if 0                /* not used */
  641. static char *plot_type_names[4] =
  642. {
  643.     "Function", "Data", "3D Function", "3d data"
  644. };
  645. static char *plot_style_names[14] =
  646. {
  647.     "Lines", "Points", "Impulses", "LinesPoints", "Dots", "XErrorbars",
  648.  "YErrorbars", "XYErrorbars", "BoxXYError", "Boxes", "Boxerror", "Steps",
  649.     "FSteps", "Vector"
  650. };
  651. static char *plot_smooth_names[5] =
  652. {
  653.     "None", "Unique", "CSplines", "ACSplines", "Bezier", "SBezier"
  654. };
  655.  
  656. static void print_points(curve)
  657. int curve;            /* which curve to print */
  658. {
  659.     register struct curve_points *this_plot;
  660.     int i;
  661.  
  662.     if (curve < 0) {
  663.     for (this_plot = first_plot, i = 0;
  664.          this_plot != NULL;
  665.          i++, this_plot = this_plot->next_cp) {
  666.         printf("Curve %d:\n", i);
  667.         if ((int) this_plot->plot_type >= 0 && (int) (this_plot->plot_type) < 4)
  668.         printf("Plot type %d: %s\n", (int) (this_plot->plot_type),
  669.                plot_type_names[(int) (this_plot->plot_type)]);
  670.         else
  671.         printf("Plot type %d: BAD\n", (int) (this_plot->plot_type));
  672.         if ((int) this_plot->plot_style >= 0 && (int) (this_plot->plot_style) < 14)
  673.         printf("Plot style %d: %s\n", (int) (this_plot->plot_style),
  674.                plot_style_names[(int) (this_plot->plot_style)]);
  675.         else
  676.         printf("Plot style %d: BAD\n", (int) (this_plot->plot_style));
  677.         if ((int) this_plot->plot_smooth >= 0 && (int) (this_plot->plot_smooth) < 6)
  678.         printf("Plot smooth style %d: %s\n", (int) (this_plot->plot_style),
  679.                plot_smooth_names[(int) (this_plot->plot_smooth)]);
  680.         else
  681.         printf("Plot smooth style %d: BAD\n", (int) (this_plot->plot_smooth));
  682.         printf("Plot title: '%s'\n", this_plot->title);
  683.         printf("Line type %d\n", this_plot->line_type);
  684.         printf("Point type %d\n", this_plot->point_type);
  685.         printf("max points %d\n", this_plot->p_max);
  686.         printf("current points %d\n", this_plot->p_count);
  687.         printf("\n");
  688.     }
  689.     } else {
  690.     for (this_plot = first_plot, i = 0;
  691.          i < curve && this_plot != NULL;
  692.          i++, this_plot = this_plot->next_cp);
  693.     if (this_plot == NULL)
  694.         printf("Curve %d does not exist; list has %d curves\n", curve, i);
  695.     else {
  696.         printf("Curve %d, %d points\n", curve, this_plot->p_count);
  697.         for (i = 0; i < this_plot->p_count; i++) {
  698.         printf("%c x=%g y=%g z=%g xlow=%g xhigh=%g ylow=%g yhigh=%g\n",
  699.                this_plot->points[i].type == INRANGE ? 'i'
  700.                : this_plot->points[i].type == OUTRANGE ? 'o'
  701.                : 'u',
  702.                this_plot->points[i].x,
  703.                this_plot->points[i].y,
  704.                this_plot->points[i].z,
  705.                this_plot->points[i].xlow,
  706.                this_plot->points[i].xhigh,
  707.                this_plot->points[i].ylow,
  708.                this_plot->points[i].yhigh);
  709.         }
  710.         printf("\n");
  711.     }
  712.     }
  713. }
  714. #endif /* not used */
  715.  
  716. static void print_table(this_plot, plot_num)
  717. struct curve_points *this_plot;
  718. int plot_num;
  719. {
  720.     int i, curve;
  721.  
  722.     for (curve = 0; curve < plot_num;
  723.      curve++, this_plot = this_plot->next_cp) {
  724.     fprintf(gpoutfile, "#Curve %d, %d points\n#x y type\n", curve, this_plot->p_count);
  725.     for (i = 0; i < this_plot->p_count; i++) {
  726.         fprintf(gpoutfile, "%g %g %c\n",
  727.             this_plot->points[i].x,
  728.             this_plot->points[i].y,
  729.             this_plot->points[i].type == INRANGE ? 'i'
  730.             : this_plot->points[i].type == OUTRANGE ? 'o'
  731.             : 'u');
  732.     }
  733.     fputc('\n', gpoutfile);
  734.     }
  735. /* two blank lines between plots in table output */
  736.     fputc('\n', gpoutfile);
  737.     fflush(gpoutfile);
  738. }
  739.  
  740. /*
  741.  * This parses the plot command after any range specifications. To support
  742.  * autoscaling on the x axis, we want any data files to define the x range,
  743.  * then to plot any functions using that range. We thus parse the input
  744.  * twice, once to pick up the data files, and again to pick up the functions.
  745.  * Definitions are processed twice, but that won't hurt.
  746.  * div - okay, it doesn't hurt, but every time an option as added for
  747.  * datafiles, code to parse it has to be added here. Change so that
  748.  * we store starting-token in the plot structure.
  749.  */
  750. static void eval_plots()
  751. {
  752.     register int i;
  753.     register struct curve_points *this_plot, **tp_ptr;
  754.  
  755.     int some_functions = 0;
  756.     int plot_num, line_num, point_num, xparam = 0;
  757.     char *xtitle;
  758.     int begin_token = c_token;    /* so we can rewind for second pass */
  759.  
  760.     int uses_axis[AXIS_ARRAY_SIZE];
  761.  
  762.     uses_axis[FIRST_X_AXIS] =
  763.     uses_axis[FIRST_Y_AXIS] =
  764.     uses_axis[SECOND_X_AXIS] =
  765.     uses_axis[SECOND_Y_AXIS] = 0;
  766.  
  767.     /* Reset first_plot. This is usually done at the end of this function.
  768.        If there is an error within this function, the memory is left allocated,
  769.        since we cannot call cp_free if the list is incomplete. Making sure that
  770.        the list structure is always vaild requires some rewriting */
  771.     first_plot = NULL;
  772.  
  773.     tp_ptr = &(first_plot);
  774.     plot_num = 0;
  775.     line_num = 0;        /* default line type */
  776.     point_num = 0;        /* default point type */
  777.  
  778.     xtitle = NULL;
  779.  
  780.     /*** First Pass: Read through data files ***
  781.      * This pass serves to set the xrange and to parse the command, as well
  782.      * as filling in every thing except the function data. That is done after
  783.      * the xrange is defined.
  784.      */
  785.     while (TRUE) {
  786.     if (END_OF_COMMAND)
  787.         int_error("function to plot expected", c_token);
  788.  
  789.     if (is_definition(c_token)) {
  790.         define();
  791.     } else {
  792.         int x_axis = 0, y_axis = 0;
  793.         int specs;
  794.  
  795.         /* for datafile plot, record datafile spec for title */
  796.         int start_token = c_token, end_token;
  797.  
  798.         plot_num++;
  799.  
  800.         if (isstring(c_token)) {    /* data file to plot */
  801.  
  802.         if (parametric && xparam)
  803.             int_error("previous parametric function not fully specified", c_token);
  804.  
  805.         if (*tp_ptr)
  806.             this_plot = *tp_ptr;
  807.         else {        /* no memory malloc()'d there yet */
  808.             this_plot = cp_alloc(MIN_CRV_POINTS);
  809.             *tp_ptr = this_plot;
  810.         }
  811.         this_plot->plot_type = DATA;
  812.         this_plot->plot_style = data_style;
  813.         this_plot->plot_smooth = NONE;
  814.  
  815.         specs = df_open(NCOL);    /* up to NCOL cols */
  816.         /* this parses data-file-specific modifiers only */
  817.         /* we'll sort points when we know style, if necessary */
  818.         if (df_binary)
  819.             int_error("2d binary files not yet supported", c_token);
  820.  
  821.         this_plot->token = end_token = c_token - 1;    /* include modifiers in default title */
  822.  
  823.         } else {
  824.  
  825.         /* function to plot */
  826.  
  827.         some_functions = 1;
  828.         if (parametric)    /* working on x parametric function */
  829.             xparam = 1 - xparam;
  830.         if (*tp_ptr) {
  831.             this_plot = *tp_ptr;
  832.             cp_extend(this_plot, samples + 1);
  833.         } else {    /* no memory malloc()'d there yet */
  834.             this_plot = cp_alloc(samples + 1);
  835.             *tp_ptr = this_plot;
  836.         }
  837.         this_plot->plot_type = FUNC;
  838.         this_plot->plot_style = func_style;
  839.         dummy_func = &plot_func;
  840.         plot_func.at = temp_at();
  841.         dummy_func = NULL;
  842.         /* ignore it for now */
  843.         end_token = c_token - 1;
  844.         }            /* end of IS THIS A FILE OR A FUNC block */
  845.  
  846.  
  847.         /*  deal with smooth */
  848.         if (almost_equals(c_token, "s$mooth")) {
  849.  
  850.         if (END_OF_COMMAND)
  851.             int_error("expecting smooth parameter", c_token);
  852.         else {
  853.             c_token++;
  854.             if (almost_equals(c_token, "u$nique"))
  855.             this_plot->plot_smooth = UNIQUE;
  856.             else if (almost_equals(c_token, "a$csplines"))
  857.             this_plot->plot_smooth = ACSPLINES;
  858.             else if (almost_equals(c_token, "c$splines"))
  859.             this_plot->plot_smooth = CSPLINES;
  860.             else if (almost_equals(c_token, "b$ezier"))
  861.             this_plot->plot_smooth = BEZIER;
  862.             else if (almost_equals(c_token, "s$bezier"))
  863.             this_plot->plot_smooth = SBEZIER;
  864.             else
  865.             int_error("expecting 'unique', 'acsplines', 'csplines', 'bezier' or 'sbezier'", c_token);
  866.         }
  867.         this_plot->plot_style = LINES;
  868.         c_token++;    /* skip format */
  869.         }
  870.         /* look for axes/axis */
  871.  
  872.         if (almost_equals(c_token, "ax$es") || almost_equals(c_token, "ax$is")) {
  873.         if (parametric && xparam)
  874.             int_error("previous parametric function not fully specified", c_token);
  875.  
  876.         if (equals(++c_token, "x1y1")) {
  877.             x_axis = FIRST_X_AXIS;
  878.             y_axis = FIRST_Y_AXIS;
  879.             ++c_token;
  880.         } else if (equals(c_token, "x2y2")) {
  881.             x_axis = SECOND_X_AXIS;
  882.             y_axis = SECOND_Y_AXIS;
  883.             ++c_token;
  884.         } else if (equals(c_token, "x1y2")) {
  885.             x_axis = FIRST_X_AXIS;
  886.             y_axis = SECOND_Y_AXIS;
  887.             ++c_token;
  888.         } else if (equals(c_token, "x2y1")) {
  889.             x_axis = SECOND_X_AXIS;
  890.             y_axis = FIRST_Y_AXIS;
  891.             ++c_token;
  892.         } else
  893.             int_error("axes must be x1y1, x1y2, x2y1 or x2y2", c_token);
  894.         } else {
  895.         x_axis = FIRST_X_AXIS;
  896.         y_axis = FIRST_Y_AXIS;
  897.         }
  898.  
  899.  
  900.         this_plot->x_axis = x_axis;
  901.         this_plot->y_axis = y_axis;
  902.  
  903.         /* we can now do some checks that we deferred earlier */
  904.  
  905.         if (this_plot->plot_type == DATA) {
  906.         if (!(uses_axis[x_axis] & 1) && autoscale_lx) {
  907.             if (auto_array[x_axis] & 1)
  908.             min_array[x_axis] = VERYLARGE;
  909.             if (auto_array[x_axis] & 2)
  910.             max_array[x_axis] = -VERYLARGE;
  911.         }
  912.         if (datatype[x_axis] == TIME) {
  913.             if (specs < 2)
  914.             int_error("Need full using spec for x time data", c_token);
  915.             df_timecol[0] = 1;
  916.         }
  917.         if (datatype[y_axis] == TIME) {
  918.             if (specs < 1)
  919.             int_error("Need using spec for y time data", c_token);
  920.             df_timecol[y_axis] = 1;    /* need other cols, but I'm lazy */
  921.         }
  922.         uses_axis[x_axis] |= 1;        /* separate record of datafile and func */
  923.         uses_axis[y_axis] |= 1;
  924.         } else if (!parametric || !xparam) {
  925.         /* for x part of a parametric function, axes are possibly wrong */
  926.         uses_axis[x_axis] |= 2;        /* separate record of data and func */
  927.         uses_axis[y_axis] |= 2;
  928.         }
  929.         if (almost_equals(c_token, "t$itle")) {
  930.         if (parametric) {
  931.             if (xparam)
  932.             int_error("\"title\" allowed only after parametric function fully specified", c_token);
  933.             else if (xtitle != NULL)
  934.             xtitle[0] = '\0';    /* Remove default title . */
  935.         }
  936.         c_token++;
  937.         if (isstring(c_token)) {
  938.             m_quote_capture(&(this_plot->title), c_token, c_token);
  939.         } else {
  940.             int_error("expecting \"title\" for plot", c_token);
  941.         }
  942.         c_token++;
  943.         } else if (almost_equals(c_token, "not$itle")) {
  944.         if (xtitle != NULL)
  945.             xtitle[0] = '\0';
  946.         c_token++;
  947.         } else {
  948.         m_capture(&(this_plot->title), start_token, end_token);
  949.         if (xparam)
  950.             xtitle = this_plot->title;
  951.         }
  952.  
  953.  
  954.         if (almost_equals(c_token, "w$ith")) {
  955.         if (parametric && xparam)
  956.             int_error("\"with\" allowed only after parametric function fully specified", c_token);
  957.         this_plot->plot_style = get_style();
  958.         }
  959.         /* pick up line/point specs
  960.          * - point spec allowed if style uses points, ie style&2 != 0
  961.          * - keywords for lt and pt are optional
  962.          */
  963.         LP_PARSE(this_plot->lp_properties, 1, this_plot->plot_style & 2,
  964.              line_num, point_num);
  965.  
  966.         /* allow old-style syntax too - ignore case lt 3 4 for example */
  967.         if (!equals(c_token, ",") && !END_OF_COMMAND) {
  968.         struct value t;
  969.         this_plot->lp_properties.l_type =
  970.             this_plot->lp_properties.p_type = (int) real(const_express(&t)) - 1;
  971.  
  972.         if (!equals(c_token, ",") && !END_OF_COMMAND)
  973.             this_plot->lp_properties.p_type = (int) real(const_express(&t)) - 1;
  974.         }
  975.         if (!xparam) {
  976.         if (this_plot->plot_style & 2)    /* style includes points */
  977.             ++point_num;
  978.         ++line_num;
  979.         }
  980.         if (this_plot->plot_type == DATA) {
  981.         /* actually get the data now */
  982.         get_data(this_plot);
  983.  
  984.         /* sort */
  985.         switch (this_plot->plot_smooth) {    /* sort and average, if */
  986.         case UNIQUE:    /* the style requires   */
  987.         case CSPLINES:
  988.         case ACSPLINES:
  989.         case SBEZIER:
  990.             sort_points(this_plot);
  991.             cp_implode(this_plot);
  992.             break;
  993.         default:
  994.             ;        /* keep gcc -Wall happy */
  995.         }
  996.         switch (this_plot->plot_smooth) {    /* create new data set     */
  997.         case SBEZIER:    /* by evaluation of        */
  998.         case BEZIER:    /* interpolation routines  */
  999.         case ACSPLINES:
  1000.         case CSPLINES:
  1001.             gen_interp(this_plot);
  1002.             break;
  1003.         default:
  1004.             ;        /* keep gcc -Wall happy */
  1005.         }
  1006.  
  1007.         /* now that we know the plot style, adjust the x- and yrange */
  1008.         /* adjust_range(this_plot); no longer needed */
  1009.         }
  1010.         this_plot->token = c_token;        /* save end of plot for second pass */
  1011.         tp_ptr = &(this_plot->next_cp);
  1012.  
  1013.     }            /* !is_defn */
  1014.  
  1015.     if (equals(c_token, ","))
  1016.         c_token++;
  1017.     else
  1018.         break;
  1019.     }
  1020.  
  1021.     if (parametric && xparam)
  1022.     int_error("parametric function not fully specified", NO_CARET);
  1023.  
  1024.  
  1025. /*** Second Pass: Evaluate the functions ***/
  1026.     /*
  1027.      * Everything is defined now, except the function data. We expect no
  1028.      * syntax errors, etc, since the above parsed it all. This makes the code
  1029.      * below simpler. If autoscale_ly, the yrange may still change.
  1030.      * we stored last token of each plot, so we dont need to do everything again
  1031.      */
  1032.  
  1033.     /* give error if xrange badly set from missing datafile error
  1034.      * parametric or polar fns can still affect x ranges
  1035.      */
  1036.  
  1037.     if (!parametric && !polar) {
  1038.     if (min_array[FIRST_X_AXIS] == VERYLARGE ||
  1039.         max_array[FIRST_X_AXIS] == -VERYLARGE)
  1040.         int_error("x range is invalid", c_token);
  1041.     /* check that xmin -> xmax is not too small */
  1042.     fixup_range(FIRST_X_AXIS, "x");
  1043.  
  1044.     if (uses_axis[SECOND_X_AXIS] & 1) {
  1045.         /* some data plots with x2 */
  1046.         if (min_array[SECOND_X_AXIS] == VERYLARGE ||
  1047.         max_array[SECOND_X_AXIS] == -VERYLARGE)
  1048.         int_error("x2 range is invalid", c_token);
  1049.         /* check that x2min -> x2max is not too small */
  1050.         fixup_range(SECOND_X_AXIS, "x2");
  1051.     } else if (auto_array[SECOND_X_AXIS]) {
  1052.         /* copy x1's range */
  1053.         if (auto_array[SECOND_X_AXIS] & 1)
  1054.         min_array[SECOND_X_AXIS] = min_array[FIRST_X_AXIS];
  1055.         if (auto_array[SECOND_X_AXIS] & 2)
  1056.         max_array[SECOND_X_AXIS] = max_~Ll¬    IS];
  1057.     }
  1058.     }
  1059.     if (some_functions) {
  1060.  
  1061.     /* call the controlled variable t, since x_min can also mean smallest x */
  1062.     double t_min, t_max, t_step;
  1063.  
  1064.     if (parametric || polar) {
  1065.         if (!(uses_axis[FIRST_X_AXIS] & 1)) {
  1066.         /* these have not yet been set to full width */
  1067.         if (auto_array[FIRST_X_AXIS] & 1)
  1068.             min_array[FIRST_X_AXIS] = VERYLARGE;
  1069.         if (auto_array[FIRST_X_AXIS] & 2)
  1070.             max_array[FIRST_X_AXIS] = -VERYLARGE;
  1071.         }
  1072.         if (!(uses_axis[SECOND_X_AXIS] & 1)) {
  1073.         if (auto_array[SECOND_X_AXIS] & 1)
  1074.             min_array[SECOND_X_AXIS] = VERYLARGE;
  1075.         if (auto_array[SECOND_X_AXIS] & 2)
  1076.             max_array[SECOND_X_AXIS] = -VERYLARGE;
  1077.         }
  1078.     }
  1079. #define SET_DUMMY_RANGE(AXIS) \
  1080. do{ assert(!polar && !parametric); \
  1081.  if (log_array[AXIS]) {\
  1082.   if (min_array[AXIS] <= 0.0 || max_array[AXIS] <= 0.0)\
  1083.    int_error("x/x2 range must be greater than 0 for log scale!", NO_CARET);\
  1084.   t_min = log(min_array[AXIS])/log_base_array[AXIS]; t_max = log(max_array[AXIS])/log_base_array[AXIS];\
  1085.  } else {\
  1086.   t_min = min_array[AXIS]; t_max = max_array[AXIS];\
  1087.  }\
  1088.  t_step = (t_max - t_min) / (samples - 1); \
  1089. }while(0)
  1090.  
  1091.     if (parametric || polar) {
  1092.         t_min = min_array[T_AXIS];
  1093.         t_max = max_array[T_AXIS];
  1094.         t_step = (t_max - t_min) / (samples - 1);
  1095.     }
  1096.     /* else we'll do it on each plot */
  1097.  
  1098.     tp_ptr = &(first_plot);
  1099.     plot_num = 0;
  1100.     this_plot = first_plot;
  1101.     c_token = begin_token;    /* start over */
  1102.  
  1103.     /* Read through functions */
  1104.     while (TRUE) {
  1105.         if (is_definition(c_token)) {
  1106.         define();
  1107.         } else {
  1108.         int x_axis = this_plot->x_axis;
  1109.         int y_axis = this_plot->y_axis;
  1110.  
  1111.         plot_num++;
  1112.         if (!isstring(c_token)) {    /* function to plot */
  1113.             if (parametric) {    /* toggle parametric axes */
  1114.             xparam = 1 - xparam;
  1115.             }
  1116.             dummy_func = &plot_func;
  1117.             plot_func.at = temp_at();    /* reparse function */
  1118.  
  1119.             if (!parametric && !polar) {
  1120.             SET_DUMMY_RANGE(x_axis);
  1121.             }
  1122.             for (i = 0; i < samples; i++) {
  1123.             double temp;
  1124.             struct value a;
  1125.             double t = t_min + i * t_step;
  1126.             /* parametric/polar => NOT a log quantity */
  1127.             double x = (!parametric && !polar &&
  1128.                     log_array[x_axis]) ? pow(base_array[x_axis], t) : t;
  1129.  
  1130.             (void) Gcomplex(&plot_func.dummy_values[0], x, 0.0);
  1131.             evaluate_at(plot_func.at, &a);
  1132.  
  1133.             if (undefined || (fabs(imag(&a)) > zero)) {
  1134.                 this_plot->points[i].type = UNDEFINED;
  1135.                 continue;
  1136.             }
  1137.             temp = real(&a);
  1138.  
  1139.             this_plot->points[i].z = -1.0;    /* width of box not specified */
  1140.             this_plot->points[i].type = INRANGE;    /* for the moment */
  1141.  
  1142.             if (parametric) {
  1143.                 /* we cannot do range-checking now, since for
  1144.                  * the x function we did not know which axes
  1145.                  * we were using
  1146.                  * DO NOT TAKE LOGS YET - do it in parametric_fixup
  1147.                  */
  1148.                 this_plot->points[i].x = t;        /* ignored, actually... */
  1149.                 this_plot->points[i].y = temp;
  1150.             } else if (polar) {
  1151.                 double y;
  1152.                 if (!(autoscale_r & 2) && temp > rmax)
  1153.                 this_plot->points[i].type = OUTRANGE;
  1154.                 if (!(autoscale_r & 1))
  1155.                 temp -= rmin;
  1156.                 y = temp * sin(x * ang2rad);
  1157.                 x = temp * cos(x * ang2rad);
  1158.                 temp = y;
  1159.                 STORE_WITH_LOG_AND_FIXUP_RANGE(this_plot->points[i].x, x, this_plot->points[i].type,
  1160.                   x_axis, NOOP, goto come_here_if_undefined);
  1161.                 STORE_WITH_LOG_AND_FIXUP_RANGE(this_plot->points[i].y, y, this_plot->points[i].type,
  1162.                   y_axis, NOOP, goto come_here_if_undefined);
  1163.             } else {    /* neither parametric or polar */
  1164.                 /* If non-para, it must be INRANGE */
  1165.                 this_plot->points[i].x = t;        /* logscale ? log(x) : x */
  1166.  
  1167.                 STORE_WITH_LOG_AND_FIXUP_RANGE(this_plot->points[i].y, temp, this_plot->points[i].type,
  1168.                                y_axis + (x_axis - y_axis) * xparam, NOOP, goto come_here_if_undefined);
  1169.  
  1170.               come_here_if_undefined:    /* could not use a continue in this case */
  1171.                 ;    /* ansi requires a statement after a label */
  1172.             }
  1173.  
  1174.             }        /* loop over samples */
  1175.             this_plot->p_count = i;    /* samples */
  1176.         }
  1177.         c_token = this_plot->token;    /* skip all modifers func / whole of data plots */
  1178.  
  1179.         tp_ptr = &(this_plot->next_cp);        /* used below */
  1180.         this_plot = this_plot->next_cp;
  1181.         }
  1182.  
  1183.         if (equals(c_token, ","))
  1184.         c_token++;
  1185.         else
  1186.         break;
  1187.     }
  1188.  
  1189.     if (parametric) {
  1190.         /* Now actually fix the plot pairs to be single plots */
  1191.         /* also fixes up polar&¶metric fn plots */
  1192.         parametric_fixup(first_plot, &plot_num);
  1193.         /* we omitted earlier check for range too small */
  1194.         fixup_range(FIRST_X_AXIS, "x");
  1195.         if (uses_axis[SECOND_X_AXIS]) {
  1196.         fixup_range(SECOND_X_AXIS, "x2");
  1197.         }
  1198.     }
  1199.     }                /* some_functions */
  1200.     /* throw out all curve_points at end of list, that we don't need  */
  1201.     cp_free(*tp_ptr);
  1202.     *tp_ptr = NULL;
  1203.  
  1204.  
  1205.     /* if first_plot is NULL, we have no functions or data at all. This can
  1206.        happen, if you type "plot x=5", since x=5 is a variable assignment */
  1207.  
  1208.     if (plot_num == 0 || first_plot == NULL) {
  1209.     int_error("no functions or data to plot", c_token);
  1210.     }
  1211.  
  1212.  
  1213.     if (uses_axis[FIRST_X_AXIS]) {
  1214.     if (max_array[FIRST_X_AXIS] == -VERYLARGE ||
  1215.         min_array[FIRST_X_AXIS] == VERYLARGE)
  1216.         int_error("all points undefined!", NO_CARET);
  1217.     FIXUP_RANGE_FOR_LOG(FIRST_X_AXIS, x);
  1218.     }
  1219.     if (uses_axis[SECOND_X_AXIS]) {
  1220.     if (max_array[SECOND_X_AXIS] == -VERYLARGE ||
  1221.         min_array[SECOND_X_AXIS] == VERYLARGE)
  1222.         int_error("all points undefined!", NO_CARET);
  1223.     FIXUP_RANGE_FOR_LOG(SECOND_X_AXIS, x2);
  1224.     } else {
  1225.     assert(uses_axis[FIRST_X_AXIS]);
  1226.     if (auto_array[SECOND_X_AXIS] & 1)
  1227.         min_array[SECOND_X_AXIS] = min_array[FIRST_X_AXIS];
  1228.     if (auto_array[SECOND_X_AXIS] & 2)
  1229.         max_array[SECOND_X_AXIS] = max_array[FIRST_X_AXIS];
  1230.     if (! auto_array[SECOND_X_AXIS])
  1231.         FIXUP_RANGE_FOR_LOG(SECOND_X_AXIS, x2);
  1232.     }
  1233.     if (!uses_axis[FIRST_X_AXIS]) {
  1234.     assert(uses_axis[SECOND_X_AXIS]);
  1235.     if (auto_array[FIRST_X_AXIS] & 1)
  1236.         min_array[FIRST_X_AXIS] = min_array[SECOND_X_AXIS];
  1237.     if (auto_array[FIRST_X_AXIS] & 2)
  1238.         max_array[FIRST_X_AXIS] = max_array[SECOND_X_AXIS];
  1239.     }
  1240.  
  1241.  
  1242.     if (uses_axis[FIRST_Y_AXIS]) {
  1243.     if (max_array[FIRST_Y_AXIS] == -VERYLARGE ||
  1244.         min_array[FIRST_Y_AXIS] == VERYLARGE)
  1245.         int_error("all points undefined!", NO_CARET);
  1246.     fixup_range(FIRST_Y_AXIS, "y");
  1247.     FIXUP_RANGE_FOR_LOG(FIRST_Y_AXIS, y);
  1248.     }
  1249.     if (uses_axis[SECOND_Y_AXIS]) {
  1250.     if (max_array[SECOND_Y_AXIS] == -VERYLARGE ||
  1251.         min_array[SECOND_Y_AXIS] == VERYLARGE)
  1252.         int_error("all points undefined!", NO_CARET);
  1253.     fixup_range(SECOND_Y_AXIS, "y2");
  1254.     FIXUP_RANGE_FOR_LOG(SECOND_Y_AXIS, y2);
  1255.     } else {
  1256.     /* else we want to copy y2 range */
  1257.     assert(uses_axis[FIRST_Y_AXIS]);
  1258.     if (auto_array[SECOND_Y_AXIS] & 1)
  1259.         min_array[SECOND_Y_AXIS] = min_array[FIRST_Y_AXIS];
  1260.     if (auto_array[SECOND_Y_AXIS] & 2)
  1261.         max_array[SECOND_Y_AXIS] = max_array[FIRST_Y_AXIS];
  1262.     /* Log() fixup is only necessary if the range was *not* copied from
  1263.      * the (already logarithmized) yrange */
  1264.     if (! auto_array[SECOND_Y_AXIS])
  1265.         FIXUP_RANGE_FOR_LOG(SECOND_Y_AXIS, y2);
  1266.     }
  1267.     if (!uses_axis[FIRST_Y_AXIS]) {
  1268.     assert(uses_axis[SECOND_Y_AXIS]);
  1269.     if (auto_array[FIRST_Y_AXIS] & 1)
  1270.         min_array[FIRST_Y_AXIS] = min_array[SECOND_Y_AXIS];
  1271.     if (auto_array[FIRST_Y_AXIS] & 2)
  1272.         max_array[FIRST_Y_AXIS] = max_array[SECOND_Y_AXIS];
  1273.     }
  1274.  
  1275. #define WRITEBACK(axis,min,max) \
  1276. if(range_flags[axis]&RANGE_WRITEBACK) \
  1277.   {if (auto_array[axis]&1) min = min_array[axis]; \
  1278.    if (auto_array[axis]&2) max = max_array[axis]; \
  1279.   }
  1280.  
  1281.     WRITEBACK(FIRST_X_AXIS, xmin, xmax)
  1282.     WRITEBACK(FIRST_Y_AXIS, ymin, ymax)
  1283.     WRITEBACK(SECOND_X_AXIS, x2min, x2max)
  1284.     WRITEBACK(SECOND_Y_AXIS, y2min, y2max)
  1285.  
  1286.     if (strcmp(term->name, "table") == 0)
  1287.     print_table(first_plot, plot_num);
  1288.     else {
  1289.     START_LEAK_CHECK();    /* check for memory leaks in this routine */
  1290.  
  1291.     /* do_plot now uses max_array[], etc */
  1292.     do_plot(first_plot, plot_num);
  1293.  
  1294.     END_LEAK_CHECK();
  1295.     }
  1296.  
  1297.     /* if we get here, all went well, so record this line for replot */
  1298.  
  1299.     if (plot_token != -1) {
  1300.     /* note that m_capture also frees the old replot_line */
  1301.     m_capture(&replot_line, plot_token, c_token - 1);
  1302.     plot_token = -1;
  1303.     }
  1304.     cp_free(first_plot);
  1305.     first_plot = NULL;
  1306. }                /* eval_plots */
  1307.  
  1308.  
  1309.  
  1310.  
  1311. static void parametric_fixup(start_plot, plot_num)
  1312. struct curve_points *start_plot;
  1313. int *plot_num;
  1314. /*
  1315.  * The hardest part of this routine is collapsing the FUNC plot types in the
  1316.  * list (which are garanteed to occur in (x,y) pairs while preserving the
  1317.  * non-FUNC type plots intact.  This means we have to work our way through
  1318.  * various lists.  Examples (hand checked): start_plot:F1->F2->NULL ==>
  1319.  * F2->NULL start_plot:F1->F2->F3->F4->F5->F6->NULL ==> F2->F4->F6->NULL
  1320.  * start_plot:F1->F2->D1->D2->F3->F4->D3->NULL ==> F2->D1->D2->F4->D3->NULL
  1321.  * 
  1322.  */
  1323. {
  1324.     struct curve_points *xp, *new_list = NULL, *free_list = NULL;
  1325.     struct curve_points **last_pointer = &new_list;
  1326.     int i, tlen, curve;
  1327.     char *new_title;
  1328.  
  1329.     /*
  1330.      * Ok, go through all the plots and move FUNC types together.  Note: this
  1331.      * originally was written to look for a NULL next pointer, but gnuplot
  1332.      * wants to be sticky in grabbing memory and the right number of items in
  1333.      * the plot list is controlled by the plot_num variable.
  1334.      * 
  1335.      * Since gnuplot wants to do this sticky business, a free_list of
  1336.      * curve_points is kept and then tagged onto the end of the plot list as
  1337.      * this seems more in the spirit of the original memory behavior than
  1338.      * simply freeing the memory.  I'm personally not convinced this sort of
  1339.      * concern is worth it since the time spent computing points seems to
  1340.      * dominate any garbage collecting that might be saved here...
  1341.      */
  1342.     new_list = xp = start_plot;
  1343.     curve = 0;
  1344.  
  1345.     while (++curve <= *plot_num) {
  1346.     if (xp->plot_type == FUNC) {
  1347.         /* Here's a FUNC parametric function defined as two parts. */
  1348.         struct curve_points *yp = xp->next_cp;
  1349.  
  1350.         --(*plot_num);
  1351.  
  1352.         assert(xp->p_count == yp->p_count);
  1353.  
  1354.         /* because syntax is   plot x(t), y(t) axes ..., only
  1355.          * the y function axes are correct
  1356.          */
  1357.  
  1358.  
  1359.         /*
  1360.          * Go through all the points assigning the y's from xp to be the x's
  1361.          * for yp. In polar mode, we need to check max's and min's as we go.
  1362.          */
  1363.  
  1364.         for (i = 0; i < yp->p_count; ++i) {
  1365.         if (polar) {
  1366.             double r = yp->points[i].y;
  1367.             double t = xp->points[i].y * ang2rad;
  1368.             double x, y;
  1369.             if (!(autoscale_r & 2) && r > rmax)
  1370.             yp->points[i].type = OUTRANGE;
  1371.             if (!(autoscale_r & 1))
  1372.             r -= rmin;    /* store internally as if plotting r(t)-rmin */
  1373.             x = r * cos(t);
  1374.             y = r * sin(t);
  1375.             /* we hadn't done logs when we stored earlier */
  1376.             STORE_WITH_LOG_AND_FIXUP_RANGE(yp->points[i].x, x, yp->points[i].type,
  1377.                          xp->x_axis, NOOP, NOOP);
  1378.             STORE_WITH_LOG_AND_FIXUP_RANGE(yp->points[i].y, y, yp->points[i].type,
  1379.                          xp->y_axis, NOOP, NOOP);
  1380.         } else {
  1381.             double x = xp->points[i].y;
  1382.             double y = yp->points[i].y;
  1383.             STORE_WITH_LOG_AND_FIXUP_RANGE(yp->points[i].x, x,
  1384.                  yp->points[i].type, yp->x_axis, NOOP, NOOP);
  1385.             STORE_WITH_LOG_AND_FIXUP_RANGE(yp->points[i].y, y,
  1386.                  yp->points[i].type, yp->y_axis, NOOP, NOOP);
  1387.         }
  1388.         }
  1389.  
  1390.         /* Ok, fix up the title to include both the xp and yp plots. */
  1391.         if (xp->title && xp->title[0] != '\0' && yp->title) {
  1392.         tlen = strlen(yp->title) + strlen(xp->title) + 3;
  1393.         new_title = gp_alloc((unsigned long) tlen, "string");
  1394.         strcpy(new_title, xp->title);
  1395.         strcat(new_title, ", ");    /* + 2 */
  1396.         strcat(new_title, yp->title);    /* + 1 = + 3 */
  1397.         free(yp->title);
  1398.         yp->title = new_title;
  1399.         }
  1400.         /* move xp to head of free list */
  1401.         xp->next_cp = free_list;
  1402.         free_list = xp;
  1403.  
  1404.         /* append yp to new_list */
  1405.         *last_pointer = yp;
  1406.         last_pointer = &(yp->next_cp);
  1407.         xp = yp->next_cp;
  1408.  
  1409.     } else {        /* data plot */
  1410.         assert(*last_pointer == xp);
  1411.         last_pointer = &(xp->next_cp);
  1412.         xp = xp->next_cp;
  1413.     }
  1414.     }                /* loop over plots */
  1415.  
  1416.     first_plot = new_list;
  1417.  
  1418.     /* Ok, stick the free list at the end of the curve_points plot list. */
  1419.     *last_pointer = free_list;
  1420. }
  1421.